Mother templates for gravitational wave chirps 



O 



Abstract. Templates used in a search for binary black holes and neutron stars in 
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', gravitational wave interferometer data will have to be computed on-line since the 

computational storage and retrieval costs for the template bank are too expensive. 
^ . The conventional dimensionless variable T = (c 3 /Gm)t, where m is the total mass 

of a binary, in the time-domain and a not-so-conventional velocity-like variable v = 
(jrGmf) 1 / 3 in the Fourier-domain, render the phasing of the waves independent of the 
total mass of the system enabling the construction of mother templates that depend 
only on the mass ratio of a black hole binary. Use of such mother templates in a 
. template bank will bring about a reduction in computational costs up to a factor of 

10 and a saving on storage by a factor of 100. 

o 
o 

PACS numbers: 04.80.Nn, 95.55.Ym, 95.75.Pq, 97.80.Af 
Gravitational waves from binary black holes during the last few seconds of their 
inspiral and merger are good candidates for a first direct observation by interferometric 



gravitational wave detectors presently under construction. The computational costs of 



(3jT)! generating and storing templates used in a search for these waves is quite expensive. 

Conventionally, one has used a two-parameter family of templates corresponding to the 
two masses m 1 and m 2 of the component stars in a binary (equivalently, one has also 
used the total mass m = mi +m 2 and the (symmetric) mass ratio 77 = m 1 m 2 /m 2 or the 
Newtonian and post-Newtonian chirp times, etc.). However, as is well-known, general 
relativity allows us to use a dimensionless time- variable T = (c 3 /Gm)t in terms of which 
all dynamical equations become independent of the total mass but depend only on the 
mass ratio 77. Therefore, relativists need only study a one-parameter family of black hole 
binaries of different mass ratios; the conclusion they draw will be applicable for binaries 
of different total masses. Of course, in doing so one has to use an appropriate re-scaling 
of all physical scales. 

In gravitational wave data analysis one has so-far not taken advantage of this 
scaling since: (1) one deals with the detector data acquired in the real time t and 
not the adimensional time T and (2) one requires for the purpose data analysis samples 
hk = h(tk), tk = kAt, that are equally spaced in an apriori chosen time-interval At, 
corresponding to a sampling rate f s = 1/ At. It is obvious that the samples Kk = h(T K ), 
Tk = KAT, equally spaced in T yield the required samples hk only for a binary whose 
total mass is an integral multiple of m = c 3 M c /G, where M c is some characteristic mass. 
In this article we show that a one-parameter family of samples hxiff), supplemented 
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with a simple linear interpolation, are sufficiently accurate for the construction of a two- 
dimensional search template bank hk(m, rj). This brings about a saving on computational 
costs by a factor of up to 10 and storage costs by a factor 100. 
In the rest of this article we work with units c = G = 1. 



1. Time-domain waveform 



Post-Newtonian (PN) theory computes the gravitational wave (GW) flux F(v) emitted, 
and the (dimensionless) relativistic binding energy E(v) of, a compact binary as 
expansions in the gauge independent velocity v in the system. For a binary consisting of 
two non-spinning black holes, of masses mi and m 2 , in quasi-circular orbit about each 
other, the flux and energy functions are both known to order v 5 (i.e. 2.5 PN) beyond 
the quadrupole approximation M, 0]: 
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Gravitational waves are dominantly emitted at twice the orbital frequency of the system. 
The GW frequency is related to the invariant velocity v via v = (vrm/cw) 1 ^ 3 - 111 the 
quasi-circular, adiabatic approximation one uses the energy balance between the flux of 
gravitational waves lost from the system and the rate of decay of the binding energy of 
the system: F(v(t)) = —m(dE/dt). 

The energy balance equation can be used to set up a pair of coupled, non-linear, 
ordinary differential equations (ODEs) to compute the orbital phase evolution <p(t) 
of the binary during the adiabatic regime: 
d^_v^ = dv F(v) 

dt m ' dt mE'(v) ' 1 J 

where we have made use of the relation between the orbital frequency and the GW 
frequency, viz, f ib(t) = (p(t)/(2 / !t) = /gwW/2 = v 3 (t) / (2-nm) and E'{v) denotes the 
■y-derivative of E(v) : E'(v) = dE/dv. The above differential equations allow us to 
compute a two-parameter family of phasing formulas ip™ ,v , % = kAt, corresponding to 
the parameters (to, rj), of the binary that helps us to construct search templates. Until 
now this has been the method of computing and storing templates. However, the cost 
of generating and storing templates in this way is quite expensive. 

In the restricted PN approximation, where the amplitude of the waveform is kept 
to the lowest PN order, while the phase is expanded to the highest PN order known, 
the GW radiation emitted by a binary at a distance r from Earth and sensed by an 
interferometric antenna, is described by the waveform 

h{t) = ^W^ v 2 {t)cosml (4) 
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where (j)(t) = 2ip(t) is the gravitational wave phase. C is a constant that takes values in 
the range [0, 1] depending on the relative orientation of the source and the antenna. It 
has an r.m.s (averaged over all orientations and wave polarisations) value of 2/5. 

Let us introduce a dimensionless time- variable T = t/m. In terms of this new 
variable the differential equations in equation (|3]) take the form 

The solutions of these two differential equations would yield a phasing formula <p(T; rj) 
that depends, unlike the phasing formula (p(t;m,r)), on only one parameter, namely 
the mass ratio rj. The advantage of this new phasing formula is that for the purpose of 
data analysis it is sufficient to create and store a finely sampled one-parameter family 
of mother templates ip v K , where T% = KAT. The phasing formula for a binary of a 
given total mass m and mass ratio rj, namely (p™' v , can be computed from the mother 
template, supplemented by a linear interpolation, without having to solve the ODEs all 
over again: 

¥T V = (1 -x)<Pk +x¥k+i, K=[k/m], x = k/m-K. (6) 

where for any real number /, [/] denotes the integer closest to, but smaller than, /. 
This way of computing templates saves on computational costs by a factor of between 
5 (for T— approximants) and 10 (for P-approximants), depending on the signal model, 
and storage costs by a factor of 100 (the number of templates in the 'm-coordinate'). 

In solving the differential equations above one must specify the initial conditions, 
namely the orbital phase and velocity at a chosen instant of time, say the instant of 
coalescence. These will, of course, depend on the total mass of the system. However, 
by choosing the instant of coalescence to be zero one can make any explicit dependence 
on the total mass vanish. 



2. Frequency-domain waveform 



We shall show in this Section that the result we found in the previous Section for the 
time-domain phasing formula holds good in the frequency domain too. 

The stationary phase approximation to the Fourier transform of the waveform in 
equation (U) is given by j|, ^, [/J 
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where ipj = 2irft±(/)(t), an overdot denotes a derivative w.r.t. t and tf is the stationary 
point of the Fourier phase, ipj(t) when / > and ip~j{t) when / < 0, namely ip^(tf) = 0. 
Indeed, tf turns out to be the instant when the GW frequency /gw is numerically equal 
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to the Fourier frequency / : / = <j>(tf) /(2tt) = /gw(*/)- One can solve for tf by inverting 
the PN expansion of fcwit)- 

On substituting for the stationary point tf, and consistently using the available 
PN expansions of the flux and energy functions one gets the usual stationary phase 
approximation (uSPA) to the inspiral signal 

r°° E'(v) 

^ f (tf) = 2-nftc - <p c + 2 / dv(v 3 - v 3 )-^ = 2nft c - <\> c + *(/; m, rj) (8) 

£ (V) 

where Vf = = {-nrnf) 1 ^ 3 . is the time of coalescence and (f>c is the phase at t = 
both of which will be set to zero. Just as in the time-domain, the frequency-domain 
phasing is also given, in the SPA, by a pair of coupled, non-linear, ODEs: 

alibi f) , , dt(vf) m(Tim) 1 ^ E'(v f ) 1 

From the computational point of view, solving the ODE's above is more efficient than 
computing the phase algebraically using equation (H). For most binaries the uSPA is 
sufficiently accurate. One, therefore, generates the signal directly in the Fourier domain 
using equation (|9|). This is a lot quicker than solving the ODEs in equation (H) and 
then Fourier transforming. As in the time-domain, one can reduce the seemingly two- 
parameter family of differential equations in equation @ to a one-parameter family 
by using v = (^m/) 1 / 3 as a Fourier variable instead of / and a dimensionless 'time 
function' T(v) = t(v)/m. In that case the resulting equations are 

M£)_ fo2T( .„ ) = o, ^M + fM = . (io) 

dv dv t (v) 

3. Explicit fast frequency-domain signal models 

We shall now give explicit formulas for the frequency- domain signal in ^-representation 
for T-approximants. As for P-approximants it is best to use the differential equations 
in equation (|i"0D. 



*!?(f) can be found by substituting in equation (|[) for the flux and energy from 
equations ([I]) and re-expanding the ratio E'(v)/F(v) and integrating term-by- 
term. The resulting Fourier phase can be conveniently expressed as ^>(f;m,rj) = 
Y,o^k(f)i~k, where ^ are functions only of the Fourier frequency and are the so- 
called (dimensionless) PN chirp times given in terms of the binary mass parameters 
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with vq = {irnifo) 1 ' 3 , /o is a fiducial frequency (e.g., the lower cutoff of the antenna 
response) and the are given, in terms of the scaled frequency v = f / f , by: 

*o=5^' * 2 = !? * s = ^' *< = ^ (12) 
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The amplitude in (^) is a simple power- law of the form f~ 7 ^ 6 H]. 

The new representation that can speed up signal generation is obtained by working 
in the Fourier domain with the post-Newtonian expansion parameter v. Introduce a 
velocity-like variable v. Indeed, v is nothing but the parameter Vf in equation (|8]). Use 
v as the Fourier variable instead of /. On substituting / = v 3 /(Trm), and setting the 
instant of coalescence t c and the phase at (fic to zero, we find 

*(/) = E0* v *~ 5 > ( 13 ) 

fc=0 

where the chirp parameters 9k are given by 
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^4 = 7^ TTTTTT^T + TT^V + TT7V • (14) 
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The Fourier transform in equation [7| can now be written as 
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fc=0 

Note that in this form the phase depends only on the mass ratio 77 and not on the total 
mass m0. One can take advantage of this feature in lowering the computational costs 
in generating templates in the following manner: One computes a set of look-up Tables, 
or mother templates, $f(v;r)) vs v for different values 77. The Fourier phase \&(/; m, 77), 
for a binary of mass m and mass ratio 77, can simply be read off from the appropriate 
look-up table: H/(/; m, 77) = ty(v = (nmf) 1 / 3 ; 77). Moreover, the cost of storing a one- 
parameter set of tables of H/(t>; 77) should be several orders of magnitude lower than the 
cost of storing a two-parameter family of Tables H/(/; m, 77). 



4. Practical Considerations 



The codes for generating mother templates have been implemented and tested for 
accuracy and are available for a free-download p|. We find that the waveforms generated 
by the interpolation method agree well with the waveforms directly computed using 
the ODEs. However, certain practical matters should be kept in mind while using 
these codes. Firstly, it is most accurate to work with either of the time- or frequency- 
domain phases, (p(T) or ^(v), respectively, rather than the corresponding waveforms. 
In other words, one must compute and store a mother phasing formula and never a 
mother waveform. Interpolating waveforms is much harder, and more inaccurate, than 
the phase of a waveform, which is a monotonic smooth function. Secondly, the initial 
velocity Vq of a mother template required in solving the differential equations for the 
phasing tp{t) must correspond to the binary of lowest total mass m min in a search: 
Vq = (vrm min /o) 1,/3 , where f is the lower frequency cutoff of a detector. Thirdly, if the 
sampling interval in real time is At then the sampling interval AT in adimensional time 
should be chosen equal to At/m max , where m max is the maximum mass of in a search, in 
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order for the interpolating formula to work best. The last two points can be summarised 
as: the smallest velocity in a mother template should correspond to the lightest binary in 
a search and the adimensional sampling rate should correspond to the heaviest binary. 
This means that if the real sampling rate is f t , then in a search for binaries of masses 
m < 1OM one must choose a sampling rate fx — 5 x lO~ 2 (f t /lkHz)(m ma , x /lOM Q ) fT0|j . 

Finally, the integration of the ODEs must be terminated at, or just before, the 
time when the system reaches the last stable orbit defined by the velocity at which the 
energy function has an extremum, that is E f (v\ so ) = 0. The right had side of the equation 
for velocity evolution in equation (|5[) diverges at this point and integration cannot be 
continued beyond v = v\ so . At certain post-Newtonian orders and approximations the 
velocity evolution close to, but well before, the last stable orbit is not monotonic; these 
correspond to the case where the flux function F(v) has a zero before the last stable 
orbit velocity. In those cases one must stop the integration before reaching dv/dT = 0, 
to avoid the code from crashing. 

5. Conclusions 

The use of an adimensional time, well-known in the literature on black hole binaries, 
in the time-domain, or a velocity-like parameter in the Fourier-domain, allows us to 
introduce a one-parameter family of mother templates which significantly reduces the 
computational costs of generating and storing binary inspiral search templates that will 
be used in the analysis of data from the up-coming ground-based gravitational wave 
interferometers. Mother templates of binary inspiral waveforms would also speed-up 
data analysis simulations that require a large number of waveforms of different masses. 

It is not clear at the present time what implications this representation will have 
for the actual choice of search templates. It is conceivable that one will be able to define 
a 'f-transform' of the data so that the problem becomes a one-parameter search for 
binary black holes with different mass ratios. However, let us not forget that such a 
scheme will only be able to detect those signals whose instant of coalescence is tc = 0; 
in other words only those that coalesce at the end of a given data segment. One will, 
therefore, have to explicitly search for binaries with every possible instant of coalescence, 
as opposed to the presently planned scheme wherein one takes advantage of the fact 
that one can search for (almost) all instants of coalescence, in one-go, by working in 
the Fourier- domain. A new scheme, therefore, is unlikely to significantly change the 
computational costs of a search, although the method may have certain advantages over 
the conventional approach. Whether or not an alternative method based on the mother 
templates introduced in this work is yet to be seen. 

Mother templates that are independent of the total mass of a binary will be available 
at all post-Newtonian orders, in the effective one-body approach || [jj and in the 
case of eccentric and spinning binaries too. Their use should make it easier for us 
to incorporate new parameters, such as eccentricity and spins, in our search templates 
to detect compact binaries in gravitational wave data. 
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